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Abstract 

Kramers- Moyal coefficients provide a simple and easily visualized method 
with which to analyze stochastic time series, particularly nonlinear ones. 
One mechanism that can affect the estimation of the coefficients is geo- 
metric projection effects. For some biologically-inspired examples, these 
effects are predicted and explored with a non-stochastic projection oper- 
ator method, and compared with direct numerical simulation of the sys- 
tems' Langevin equations. General features and characteristics are identi- 
fied, and the utility of the Kramers-Moyal method discussed. Projections 
of a system are in general non-Markovian, but here the Kramers-Moyal 
method remains useful, and in any case the primary examples considered 
are found to be close to Markovian. 

1 Introduction 

In deterministic nonlinear systems, embedding theorems provide one method 
of reconstructing the system's phase space and thereby obtaining a qualitative 
understanding, at least, of the system's structure. In stochastic systems no 
such embedding theorems exist [I] . One approach for analyzing stochastic time 
series has been pioneered by Friedrich et al. [H |3] and continues to be devel- 
oped by Kleinhans [H [51 [5] and others. The Kramers-Moyal coefficients (in one 
dimension) , 

D^^\xo,t) ^ \im —{[x{t + t) - xor\x{t) = xo), (1) 

are a critical part of the theory of Fokker-Planck equations, but only recently 
did Friedrich et al. suggest they could be computed directly from data. The 
first and second coefficients are known as the drift and diffusion coefficients; 
the higher-order coefficients are often zero |7]. The angle brackets denote an 
ensemble average, but if the Kramers-Moyal coefficient is constant in time this 
can be replaced with a time average. If the system can be modeled by a Langevin 
equation, or, equivalently, a Fokker-Planck equation, then the drift and diffusion 
coefficients completely characterize of the system. 



Such an approach, aheady apphed to fields including neuroscience [3], car- 
diology [8] , traffic engineering [9] , finance [lOl |TT] and turbulence [2] , appears 
well suited to molecular biology, where biological structures commonly move in 
an overdamped, strongly Brownian environment [12] . Due to structural con- 
straints, however, biological systems are often constrained in their motion to 
diS^use on a manifold which is no longer Cartesian. Experimentally, meanwhile, 
position measurements are usually taken in a Cartesian framework. Due to these 
geometrical effects the Kramers-Moyal coefficients derived from measurements 
are not the same as on the particle's manifold. We show how these projected 
drift and diffusion coefficients can be used to make inferences about the original 
system's dynamics. 

For illustration, we consider three specific examples. Motivated by the diffu- 
sional search phase of myosin- V [T3l[T4] and other molecular motors, we consider 
a particle undergoing tethered diffusion with a rigid tether, so that the particle 
is constrained to a sphere. We assume its motion is recorded only in one linear 
dimension, with an experiment like the traveling wave tracking method of Cap- 
pello et al. (TS] This involves the projection of a two-dimensional system onto 
one dimension. In this example, the diffusion on the sphere is free; in a second 
example, we make it harmonic in angle from a preferred orientation ('biased 
diffusion on a sphere'). Lastly, we consider the more complicated example of 
a particle undergoing free diffusion on a sphere about a point itself undergoing 
biased diffusion on a sphere. We refer to this as 'compound diffusion'. 

We begin in section II with a simple example to illustrate the projection 
concept. After presenting a general framework for predicting geometrical effects 
in section III, we continue with the less trivial examples in sections IV and V. 
The effects of a non-Markovian projection are discussed in section VI, followed 
by conclusions in section VII. 



2 Free diffusion on a sphere 

Consider a particle diffusing on a sphere, of radius r, under potential U and 
constant local diffusion coefficient D. Raible and Engel found that the 
stochastic differential equations (SDKs, or Langevin equations) describing this 
motion in spherical polar co-ordinates {6, (j)) are 



de = {-dgU + D cot e)dt + V2Ddwe 
^"^^ dt + ^ „ dwj,. 



(2) 



sin^ d sm 



dwe and dw^ are Wiener processes and the SDE is to be interpreted under the 
Ito convention. Setting [7 = for free diffusion and projecting onto z ~ rcosO, 
we have by the Ito rule for change of variables 



dz = -2Dzdt - v/2L>(r2 - z^)dwg. (3) 
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This corresponds to drift and diffusion coefficients 



D^^^z) = -2Dz 
D^^\z)^D{r^~z^). 



(4) 



Therefore even though in local rectangular co-ordinates on the surface of the 
sphere the drift is zero and the diffusion constant ('free diffusion'), in this pro- 
jection the drift is linear, and the diffusion is an inverted (convex) parabola. 

The drift coefficient is analogous to a position-dependent force, and the dif- 
fusion to a position-dependent temperature or 'noise' term [9]. Although the 
original motion on the surface was force-free with constant diffusion coefficient, 
we see the geometric projection introduces an effective 'restoring force', propor- 
tional to — z. ft also introduces an effective noise which is maximal at z = but 
decreases to zero at z = ±r: this prevents the projected particle from moving 
past the edge of the sphere. 

The equilibrium probability distribution of this process is [7] 



This illustrates the power of the Kramers-Moyal approach. A naive analysis 
of a stationary stochastic time series would be to construct its probability dis- 
tribution. For the example above, this would give a uniform distribution, from 
which a reasonable guess for the system would be free diffusion in a one dimen- 
sional infinite potential well (corresponding to constant diffusion, and drift zero 
and infinite inside and at the edges of the well, respectively). On computing 
the drift and diffusion coefficients from a time series by Eq. ([1]), however, one 
would recover the drift and diffusion coefficients of Eqs. ([4]), and obtain a dif- 
ferent characterization. One explanation of these drift and diffusion coefficients 
could be to recognize the system as a projection of free diffusion on a sphere, 
although this identification cannot be unique since we are projecting onto one 
dimension a two-dimensional system, and without a priori knowledge could be 
three- or higher-dimensional. Finally, if the projection is Markovian, as dis- 
cussed below, the drift and diffusion coefficients fully characterize the dynamics 
in this co-ordinate, which the probability distribution does not. 

Unlike stochastic autoregressive-type models [P, we see that the Kramers- 
Moyal approach permits arbitrary nonlinearity in the drift and diffusion coeffi- 
cients to be recovered, and easily visualized. 

3 General analysis of projections 

Suppose we project the dynamics of a system x{t) onto some (collective) variable 
X{t) = f{x{t)). Suppose the system can be written as a set of first-order 




For free diffusion on a sphere, this gives a uniform distribution 




\z\ < r 
\z\ > r 



(5) 
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stochastic differential equations. With the use of Ito's chain rule, an SDE for 
X in terms of x may be derived, 



dX = gx {x, t)dt + hx {x, t)dw. (6) 

We consider the estimation of the Kramers-Moyal coefficients ([T]) for this pro- 
cess. Writing the conditional expectation of Eq. ([T|) in an alternative form 
gives 

D^ixo^t) = lim ^ / {X-Xo)'^P{X,t + T\xQ,t)dX, 
T^o n'.T J 

where Xq — /(xq). We find Dx\x,t) = gx{x,t), Dx\x,t) = ^hx{x,t), and 

D^'' (x, t) = for n > 3 [7] . But these coefficients still involve the full system 
variables x, to which, we assume, the experiment does not have access. 

From basic laws of conditional probability, it can be shown that the projec- 
tion of these Kramers-Moyal coefficients from x onto the subspace X (which we 
distinguish from the unprojected coefficients by its arguments) is 



n(")^v ^-\ [ r,(-"^)(^ 5{Xq - /(fo))P(xo, t) 
Dy{Xo,t)= / D'x'[xQ,t) p^^^ dxQ. (7) 



This constitutes a projection operator for the Kramers-Moyal coefficients, 
Vx 



D^-\x,t)\ ^ j DPix,t)'-^^i^l^^dx. 

Clearly Vx'Px = Vx so Vx is indeed a projection. As Z?(")(x, t) is a dynamical 
variable on the full phase space x, the Vx above coincides with the projection 
operator from the standard Zwanzig-Mori projection techniques }17j . To reiter- 
ate, however: before applying the projection operator to an arbitrary Langevin 
equation, one must transform the equation into the form of Q , one in the vari- 
able into which the system is to be projected. The quotient 5{X — f{x)) p\^x% 
is also known as the conditional distribution for the subset of the full phase 
space distribution P{x, t) on the surface f{x) = X. The projection also has 
implications for the Markov property, as will be discussed in section |6l 

If the 1-point densities P{xQ,t) and P{Xo, t) are well-defined, so are the pro- 
jections d'"-^\x, t) (although they in some cases may be trivially zero). Further- 
more, since Z?^^(x, <) — for n > 3, then so do D^^\x,t). The projected drift 
and diffusion coefficients, the equations above show, may be time-dependent 
even if the original system is homogeneous (the coefficients in its equations 
of motion independent of time) if the distribution P{xQ,t) is time-dependent 
(non-stationary). If the original system is both homogeneous and stationary, 
the projected drift and diffusion coefficients will be stationary. 

The Kramers-Moyal coefficients may of course be computed directly from 
experimental data ([T]) of the projected variable. The necessary averages may 
be computed by time averaging if the projected drift and diffusion coefficients 
are stationary, which as just noted is when the full system is homogeneous and 
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stationary. (This is in contrast to Kramers-Moyal analysis on the full system, 
where the system need only be homogeneous.) This projection operator method, 
however, provides a means of predicting the projected drift and diffusion coef- 
ficients from a model of the full system without the need for direct numerical 
simulation. 

We now illustrate these observations with two examples. 



4 Biased diffusion on a sphere 

Suppose a particle is again undergoing diffusion on a sphere with local diffusion 
coefficient D, but with a potential harmonic with angle away from a preferred 
binding direction. Align the z-axis with this preferred direction so that U = 
kO"^ /2. Substituting into the Langevin equations ^ for diffusion on a sphere, 
the equilibrium solution is 



Pe,j,(B) oc sinl 



normalized such that JJ P0ci,{d)dddip — 1. 

Let us project this motion onto an axis s tilted at angle 9q from the 2;-axis 
(9 — direction) in the x ^ z plane, 

s{9, (j)) — z cos 9q + X sin 9q 

= r cos^o cos6' + sin^o sin6'cos0. (8) 

The differential 

ds = gs{9, <i))dt + hs{9, (j))dw (9) 

can be calculated using Ito's law for change of co-ordinates, and adding the 
noise terms, which we assume to be independent, in quadrature. 

For the steady-state probability distribution as a function of s, we calculate 



Psiso) - ((5(so-s(0,0))) 

IS 



the second equality by integrating over (j>, and where we denote 3^(0, 0o) = 
d^s{9, The single integral that results we calculate numerically. 

Similarly, from Eq. ([7]) we project the Kramers-Moyal coefhcients by 

where we have performed the integral over (j) analytically, and in which (j)Q is 
itself a function of 9. 

In Fig. [T] we show drift and diffusion coefficients for s(t) computed from 
direct numerical simulation of Eqs. ^ and by the semi-analytical approach 
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culminating in Eq. (|10p . The agreement between numerical and theoretical 
approaches is excellent, with two exceptions. At x = ibcos^o there are spikes 
due to the 6*0 = 0,7r singularities in Eqs. Raible and Engel [H] transform 
their equations into three-dimensional Cartesian co-ordinates before performing 
computations; we retain the spherical polar co-ordinates for transparency. Sec- 
ondly, the larger binding potential for fc = 5 means the large negative regions 
of X are poorly sampled and also may lead to finite sampling time effects |18) . 

From Fig. [1] we see the potential on the sphere U — k0^/2 has an effect 
additional to the geometric effective force arising from the projection, and it is 
largest away from the potential minimum z — cos^o- The projected diffusion 
coefficients, meanwhile, show no change with the strength of the potential on 
the surface, k. We infer that only the geometry of the system has affected 
the diffusion coefficient, and therefore that an inverted parabola shape for the 
diffusion coefficient is a good indicator of diffusion on a sphere. 

More generally, if a particle is diffusing on a surface that is radially symmetric 
about the measurement axis s, with radius r(s), then the projected diffusion 
coefficient is independent of the potential on the surface. It can be shown from 
Raible and Engel's formulas that the projected diffusion coefficient will be 



assuming that the diffusion coefRcient on the surface has constant value D. 

5 Compound diffusion on a sphere 

We extend the previous scenario to one where there is a second rod undergoing 
tethered diffusion about the end of the first. We assume the first rod has a 
preferred binding angle, but that the second rod is diffusing freely. This models 
the motion of the unbound head during the diffusional search phase of molecular 
motors such as myosin-V. 

In a first approximation, we assume that the motion of the end of the first 
rod with respect to its tethered end is given by the biased diffusion on a sphere 
of section IV, and the motion of the free end of the second rod with respect 
to the end joined to the first rod is given by the free diffusion on a sphere of 
section II. Once again we project onto a linear co-ordinate. The position of the 
end of the first rod is then s{9, (/)) from Eq. ([5]), and the end of the second rod 
X{9,4',X2) = s{0,(l)) + X2, with X2 given by Eq. ^ for z. The equilibrium 
distribution of X2, Px2{x2) is therefore given by Eq. (O. We further assume 
that the processes s{9,(f)) and X2 are independent. Therefore the (unprojected) 
drift and diffusion coefficients are 




D 



1 + (dr/ds) 



2 ' 



gx{0,(t>,x2) = gs{0,4') +gx2{x2) 
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Position, X 



Figure 1: Projected (top) drift and (bottom) diffusion coefficients for biased 
diffusion on a sphere. Semi-analyticai predictions (|f Op are shown with soiid 
fines, together with resufts from direct numerical simulation ^ for restoring 
force strength fc = (circles), 1 (triangles) and 5 (squares). Preferred binding 
angfes (and reference for co-ordinate transformation) were Oq = 7r/2, 7r/4, 7r/4 
respectively. The diffusion coefficient on the surface was D — 1 and the radius 
of the sphere (iength of the tether) r = 1 in all cases. 
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Wc first require the equilibrium probability distribution for X. Integrating 
{5{X — s{9, (f>) — X2))e,^,x2 over X2 leaves the convolution 

Px{X) = JJ Pe4,{0)Px2{X - s{e,4>))ded4>, 

which can be computed numerically. For the drift and diffusion coefficients, 
from Eq. ([7]) we obtain, again after integrating over x^^ 

Pe4>{0A)Px2{X - s{e,^))d9d(l). 

Numerical results again together with semi-analytical predictions are shown 
in Fig. [21 Observe the unusual bimodal shape for the diffusion coefficient 
at A: = 0, reflecting the two 'spheres' making up the diffusion process. As k 
increases this shape begins changing towards a single inverted parabola which 
one would expect in the limit k ^ oo where the first rod is fixed. The drift 
is linear at fc = 0. It maintains an approximately linear region but at larger k 
develops large swings beyond x — cosOq ± 1 « —0.3,-1-1.7, regions which will 
eventually become forbidden with a fixed first rod (fc — s- oo). 

When the motor is in the diffusional search state, 1-D Kramers-Moyal analy- 
sis on experiments tracking the unbound head of myosin-V, or similar molecular 
motors, should (in the absence of other effects) give results like the above. Al- 
ternatively, if the marker is connected to a point on the unbound neck domain, 
the above simulations should be repeated with a shorter length for the second 
rod. Yildiz et al. [IS] have performed such experiments on kinesin, from which 
they identified the step size and type (hand over hand), but their large sampling 
interval yields insufficient data for the present analysis. 

6 Markov property 

In our examples, the full system is Markov, that is, its next state depends only 
on the current one and not any other history: 

P{x2,t2\xi,ti) = P{x2,t2\xi,ti;xo,to), 

where t2 > ti > to- Even if the full system is Markov, however, in general a 
projection is not [T7]. Then the projected variable cannot be fully modeled by a 
Langevin equation (with history-independent coefficients), and so the drift and 
diffusion coefficients, or even all the Kramers-Moyal coefficients, cannot fully 
characterize the dynamics of the projected variables. Some authors [[Tj have 
criticized the Kramers-Moyal method on this basis, since experimentally one 
invariably can only access projections, not the entire system. 

First, we wish to emphasize that even for projections that are demonstrably 
not Markovian, the Kramers-Moyal coefficients remain well-defined. In this 
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Position, X 

Figure 2: Semi-analytical and numerical estimates for projected (top) drift and 
(bottom) diffusion coefficients for compound diffusion on a sphere, as described 
in the text. Parameters for the first rod are the same as for Fig. [l] and with 
the same marker shapes. Surface diffusion coefficient and radius for the second 
rod (with respect to the first) are also D — 1 and r = 1, respectively. 
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case one may not write down an equation of motion of the system from the 
reconstructed coefficients, but they remain a useful tool for characterization 
and comparison of time series. Conversely, the existence of Kramers-Moyal 
coefficients, or 2-point probability densities, of course do not by themselves 
infer anything about the Markov property of the underlying process. 

As an example of the Kramers-Moyal coefficients in a non-Markovian pro- 
jection, consider the two-dimensional Ornstein-Uhlenbeck process 



projected onto X — xi. With a specific set of parameters, we plot the drift and 
diffusion coefficients in Fig. [3l As we will show below, this projection is non- 
Markovian. For this system the projection operation reduces to / d'"^\x, X2)P{X, X2)dx2/ P{X). 
Qualitatively, the projected drift and diffusion coefficients for a value of the pro- 
jected variable X = xi are the values of the unprojected coefficients averaged by 
the likely locations of X2 given this xi. In other words, they are the expectations 
over the conditional probability density P{x2\xi). For this system we see in Fig. 
[3] an effective deterministic force tending to restore X towards the origin, that 
is confine it to a finite region, and that there is a noise source operating directly 
on the X variable. 

There are some cases where we know the projection will be Markovian. If 
the projected variable is slow compared to the other degrees of freedom of the 
system, the projection will be Markovian [17| . In the unbiased diffusion on a 
sphere of section [5] we saw that because of the special symmetry of the system, 
the second of the two degrees of freedom {(p) was irrelevant for the projected 
variable, which by Eq. ([3|) we see is Markovian. 

The Chapman-Kolmogorov equation [7] 



although not sufficient for the Markov property, is generally taken as a good test 
of it pUl nn . We tested the biased diffusion on a sphere of section |3] with 
the Chapman-Kolmogorov equation and the a posteriori test of Michcletti et al. 
[2T| . Stationarity was assumed and the transition probabilities were calculated 
by time average, setting t2 — ti — ti — to — r. The results are shown graphically 
in Figs. IHandEl These tests indicate that this projection of biased diffusion on 
a sphere is at least approximately Markovian. To a good approximation, then, 
the projected variable can be modeled by a Langevin equation in that variable, 
and in which the noise is Gaussian and delta-correlated. 

For comparison we show in Fig. [6] clear failure of the Chapman-Kolmogorov 
equation and Micheletti tests for a non-Markovian system, that of Eq. (|lip 
projected onto X = xi. 



dxi = {axi + bx2)dt + ^2DidWi 



(11) 



dx2 = [cxi + dx2)dt + ^/2D2dW2 
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Figure 3: Sample trajectory (top) and projected (left) drift and (right) diffusion 
coefficients for the projection X — xi oi the system of Eqs. pT|) with a = 0, 
b = 1, c — —1, d = —0.1, Di = 1 and D2 = 10 by direct numerical simulation 
(circles) and analytically (lines). Analytical predictions were calculated using 
Erickson's results [19] for multi-dimensional Ornstcin-Uhlenbeck processes, giv- 
ing D(i)(X) = and D^^^X) = D^, under smtable stability 

conditions on the SDEs (fTTj) . 




Figure 4: Chapman-Kolmogorov equation for projected biased diffusion on a 
sphere Q with r = 0.1. Contours, in logarithmic scale, of the left-hand side of 
the equation are in gray and right-hand side in black. 
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Figure 5: (left) Normalized autocorrelation and (right) probability density (cir- 
cles) with Gaussian fit (solid line) of the reconstructed noise, as per the the 
a posteriori test of Micheletti et al. 21\, for projected biased diffusion on 
a sphere The test reconstructs the effective noise dW from the recon- 

structed drift and diffusion coefficients, assuming a Langevin equation of the 
form dX = g{X)dt + h{X)dW holds. If dW has delta-autocorrelation and 
Gaussian probability density, claims the test, the process is Markov. 




Figure 6: (left) Chapman-Kolmogorov equation with t = 1 for the projection 
onto X = xi the Ornstein-Uhlenbeck process pT|) with parameters as in Fig. 
[3l Contours, in logarithmic scale, of the left-hand side are in gray and right- 
hand side in black. They clearly do not coincide, (right) Autocorrelation of the 
reconstructed noise, as per the Micheletti test described in Fig. [5l The noise is 
clearly not delta-autocorrelated. 
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7 Conclusions 



We began with a simple example of projecting free diffusion on a sphere onto 
a rcctangTilar co-ordinate. A general projector operator formalism for pro- 
jecting Kramers-Moyal cocfBcicnts was then established. This was applied to 
biologically-inspired examples of biased diffusion on a sphere, and compound 
diffusion on a sphere, where excellent agreement with numerical results was 
observed. Generally speaking, the diffusion coefficient is more useful for iden- 
tifying for the geometry of the surface on which the object is diffusing, the the 
drift coefficient for the potential on the surface. It was noted that the Kramers- 
Moyal coefficients can provide useful information even when the projection is 
not Markovian, which in general is the case, although the examples considered 
here were close to Markovian. 
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